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I. Introduction 


A major difficulty when solving the incompressible flow equations that use primi- 
tive variables is caused by the pressure term which is used as a mapping parameter 
to obtain a divergence-free velocity field. One commonly used approach is to solve 
the Poisson equation for pressure, which is derived from the momentum equations 
[1]. This approach can be very time consuming. To accelerate the pressure- 
field solution and alleviate the drawback associated with the Poisson equation ap- 
proach, Chorin [2] proposed the use of artificial compressibility in solving the con- 
tinuity equation. A similar method was adopted by Steger and Kutler [3] and 
Chakravarthy [4] using an implicit approximate-factorization scheme [5]. Based on 
this procedure, a pseudocompressible method has been developed for solving three- 
dimensional, viscous, incompressible flow problems cast in generalized curvilinear coor- 
dinates [6,7]. The purpose of the present paper is to show salient features of the pseudo- 
compressible approach, which is primarily designed for obtaining steady-state solutions 
efficiently. 

' n. Description of the Method 


In the present formulation, the three-dimensional, incompressible Navier-Stokes equa- 
tions are modified to form the following set of governing equations written in dimensionless 
form : 


1 dp dUi 
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dm ^ 

dt dxj dxi dxj 


(16) 


Here, t is time; Xi are the Cartesian coordinates; u, are corresponding velocity com-, 
ponents; p is the pressure; and r,y is the viscous stress tensor. The parameter 

IjP is the pseudocompressibility. As the solution converges to a steady state, 
the pseudocompressibility effect approaches zero, yielding the incompressible form 
of the equations. In the present study, the approximate factorization scheme 
by Beam and Warming [5] is implemented to solve the finite-difference form of 
the governing equations written in general curvilinear coordinates (see ref. 6 for 
detail). 

In the present formulation, waves of finite speed are introduced. And the system 
of modified equations given by equations (la) and (lb) can be marched in time. The 
magnitude of the wave speed depends on To recover the incompressible phenomena, 



the physics requires that the pressure wave propagates much faster than the spreading 
of vorticity. From this, the following criterion for the lower bound on is obtained 

[TJ: 

;? > (1 + 4{Xref 1X6 fi^L href )IRtf “1 ( 2 ) 

where ir«/ is the reference length, and xs and xt are the characteristic lengths that the 
vorticity and the pressure waves have to propagate during a given timespan. 

The upper bound on P depends upon the particular numerical algorithm chosen. 
In the present study, higher-order cross-differencing terms are added to obtain the 
approximately factored form of the governing equations. These added terms con- 
taminate the momentum equations as well as the continuity equation, and there- 
fore must be kept smaller than the original terms everywhere in the computational 
domain. This requirement leads to the following criterion for the upper bound of 

P- 

/?At<0(1) (3) 

where Ar is the time-step used in the integration scheme. 

in. Computed Results 


Numerical experiments were performed to illustrate the present procedure. To 
represent an internal flow, the flow through a channel at Re=l,000 was chosen. 
The coordinate system and velocity vectors for a converged solution are shown in 
flgures la and lb. To change the ratio of the time scales required for the pressure 
waves and the vorticity to map the entire flow fleld, the channel length, L, is varied 
form 20 to 40. The recommended values of P for these cases using Ar = 0.1 
are: 


0.75 < Pl^2Q < 10, 1.19 < Pl=zo < 10, 1.69 < /?l=40 < 10 

In table 1, the number of iterations for one roundtrip by the pressure wave (denoted 
by Ni ) is tabulated for various values of P which include values outside the recommended 
range. In flgure 2, root-mean-square (RMS) values of {divu) are plotted to check the 
accuracy of the converged solutions. When the value of P is out of the range specifled, the 
accuracy of the solution deteriorates. 

To represent an external flow, the flow past a circular cylinder at a Re = 40 was 
chosen. To obtain the near-fleld solution only, the distance traveled by the waves and 
the spreading of the vorticity can be approximately the same in magnitude. In the 
present case, this leads to the range for P using Ar = 0.1 to be 0.1 < < 10 . 

This indicates that the magnitude of P is less restrictive for external flows. In flgures 
3a and 3b, the stream-function contours and the pressure coefficient on the surface 
are shown for a steady-state solution. This solution agrees very well with that of 
Mehta who used a stream function and vorticity formulation in two dimensions (private 
communication, U. B. Mehta, 1983). In flgure 4, in which the history of the pressure 
drag is shown for an impulsively started circular cylinder at Re = 40, four different 
values of p were compared with the time-accurate solution of Mehta. In all cases, the 
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values of P are selected within the suggested range above, and the solutions converge 
rapidly. 

To test internal flows further, an annular duct with a 180® bend is chosen. This 
conflguration is similar to the turnaround duct of the hot-gas manifold in the Space 
Shuttle main engine (SSME). In flgures 5a and 5b, the geometry and a laminar solution 
at Re= 1,000 are shown, which reveals the formation of a large separated bubble after 
the 180® bend. For this geometry, the streamwise length normalized by the duct 
width is 20. The test problems presented here were treated using a 51 x 17 x 21 
mesh for half-duct formulation and the computing time required was 1.1 x 10“^ sec 
per mesh point per time-step on the Cray X-MP computer at NASA Ames Research 
Center. 


IV. Concluding Remarks 

This paper presents salient features of the computational procedure developed for a 
three-dimensional, incompressible, Navier-Stokes code. This procedure has been applied 
to various geometrically complex flows, including a major application in analyzing the 
flow fleld in the SSME power head. The present algorithm has been shown to be very 
robust and accurate if the selection of is made according to the guidelines presented 
here. 
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(b) Fully developed velocity profile 


Figure 1.- Developing laminar channel flow at 
Re = 1,000 (Re based on channel width and 
average velocity). 
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Figure 2.- RMS (divu) history of chsimel flow at 
Re = 1,000 and At = 0.1. 
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